Smooth Phases, Roughening Transitions and Novel Exponents 
in One-dimensional Growth Models. 
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A class of solid-on-solid growth models with short range interactions and sequential updates is stud- 
ied. The models exhibit both smooth and rough phases in dimension d = 1. Some of the features 
of the roughening transition which takes place in these models are related to contact processes or 
directed percolation type problems. The models are analyzed using a mean field approximation, 
scaling arguments and numerical simulations. In the smooth phase the symmetry of the underlying 
dynamics is spontaneously broken. A family of order parameters which are not conserved by the 
dynamics is defined as well as conjugate fields which couple to these order parameters. The corre- 
sponding critical behavior is studied and novel exponents identified and measured. We also show 
how continuous symmetries can be broken in one dimension. A field theory appropriate for studying 
the roughening transition is introduced and discussed. 
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I. INTRODUCTION 

The statistical properties of moving interfaces and sur- 
faces of growing crystals have been studied extensively 
in recent years 0. Various theoretical approaches have 
been applied in these studies including field theoreti- 
cal analyses of continuum KPZ-type equations J|-§] and 
studies of discrete growth processes such as solid-on-solid 
(SOS) or polynuclear growth (PNG) models amongst 
others |-@. 

One of the more interesting properties of moving in- 
terfaces is their roughness. It is well known that in d > 2 
dimensions, moving interfaces may be either rough or 
smooth depending on the level of the noise in the system. 
However, growth processes in d = 1 dimensions are more 
subtle. Most growth processes governed by short range 
interactions, such as those described by the KPZ equa- 
tion yield a rough interface. On the other hand study of 
a class of PNG models suggested that both smooth and 
rough phases may exist in one dimension (Id) |p|, p^| , |l6| . 
However, this class of models is characterized by two key 
features: (a) The evolution takes place by a parallel up- 
dating process in which time progresses in discrete steps 
and all sites of the interface are updated simultaneously 
according to a given rule at any given time step. Such 
dynamics tend to be less noisy than sequential updat- 
ing processes, in which one site is updated at a time, 
(b) The models have a maximal velocity by which the 
uppermost point of the surface can propagate. The exis- 
tence of a maximal velocity in these models is due to the 
use of parallel updates, and the smooth phase disappears 
when random sequential (continuous time) updates are 
used. An interesting question is whether a sequential up- 
date growth process is capable of exhibiting smooth and 
rough phases in Id. 



Some time ago a transition from a smooth to a rough 
surface was observed in a SOS model for surface recon- 
struction with sequential updates and particle conserva- 
tion [jrH| . The free parameter, if, is the maximal allowed 
height difference between adjacent sites and is a discrete 
quantity. It was observed that the surface is smooth for 
H < 3 and macroscopically grooved for if > 3. At if = 3 
the surface appears to be rough. This phenomenon is 
clearly related to the local conservation law in this model. 

Recently, a class of growth processes with short range 
interactions, sequential updates and non-conserved dy- 
namics was introduced |jl7). They were demonstrated 
to exhibit both rough and smooth phases in Id. The 
dynamics is described by SOS models in which adsorp- 
tion and desorption processes take place in a ring geome- 
try. Depending on the relative rates of the two processes 
one finds either a smooth or a rough phase. In studying 
the roughening transition in these models, it has been 
found that some of its features are related to those of 
contact processes, or directed percolation, which have 
been extensively studied in the past |18|-pi|. These lat- 
ter systems exhibit a continuous phase transition which 
is strongly linked to the existence of absorbing states (a 
set of states from which the system can not escape) . The 
model introduced in (Tij ] may be viewed as composed of 
a series of contact processes interacting with each other, 
whose dynamics does not have an absorbing state. This 
model is the subject of the present paper. 

An intriguing question related to the existence of a 
smooth phase is that of spontaneous symmetry breaking 
(SSB) and long-range order in Id systems with short- 
range interactions and small but unbounded noise. By 
unbounded noise we mean that in a finite system all mi- 
croscopic configurations can be reached from any initial 
condition in a finite time. It is well known that in thcr- 
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mal equilibrium no phase transition takes place under 
these conditions. Systems far from equilibrium p2| such 
as moving interfaces or driven diffusive systems are, how- 
ever, less restrictive and the question of existence of SSB 
under these conditions has been open for some time. Re- 
cently, a simple non-equilibrium model which exhibits 
SSB in Id was introduced p3| , p4|] . This model belongs 
to a class of driven-diffusive systems, in which charges of 
two kinds are injected at the ends of a Id lattice and are 
biased to move in opposite directions. The dynamic rules 
are symmetric under space and charge inversion (PC in- 
variance). However, this symmetry is broken in the sta- 
tionary state of the system where the currents of the two 
charge species are different. In this model, SSB is the re- 
sult of the conserved order parameter in the bulk (charges 
are not created or annihilated except at the boundaries) , 
and the existence of open boundary conditions. These 
two conditions favor the emergence of SSB; the conserved 
dynamics slows down the temporal evolution of the sys- 
tem, moreover flips from one broken symmetry phase 
to another can originate only at the boundaries where 
the order parameter is not conserved. Initial attempts 
to modify this model such that either or both of these 
features are eliminated resulted in symmetric stationary 
states with no SSB. 

The growth model discussed in the present work pro- 
vides a simple example for SSB which takes place far 
from equilibrium. The breaking of symmetry takes place 
in the smooth phase and the order parameter associated 
with this transition is not conserved by the dynamics. 
The model thus demonstrates that SSB can take place in 
Id in non-equilibrium systems with non-conserved order 
parameter and ring geometry. 

In this paper we present a detailed analysis of the 
growth models introduced in jl7| ]. The relation of the 
models to contact processes and directed percolation is 
discussed and several families of novel critical exponents 
are defined. One family describes the critical behavior of 
magnetization-like order parameters related to the sym- 
metry breaking. Another family is related to statistical 
properties of the interface height near the roughening 
transition. The dynamical equations are analyzed us- 
ing the mean field approximation, and a field theoretical 
model appropriate for the roughening transition is in- 
troduced. The scaling picture that emerges is far from 
complete but points to the existence of complex critical 
behavior. 

The paper is organized as follows: The growth model is 
defined in Sec. 0. The relation to contact processes and 



directed percolation is discussed in Sec. III. In Sec. IV 



the results of scaling analysis and numerical studies are 
presented. The question of spontaneous symmetry break- 
ing, and the family of order parameters and their crit- 
ical behavior are discussed in Sec. [v|. The dynamics 
of the model is studied within the mean field approxi- 
mation in Sec. VI, and a field theoretical model appro- 
priate for studying the critical behavior of the roughen- 
ing transition is defined and discussed in Sec. VII. In 



Sec. VIII some light is shed on the relation of the model 
to the polynuclear growth models discussed above and in 
Sec. [lX| the main results are summarized and conclusions 
drawn. Finally, a particular case for which steady state 
can be calculated exactly is presented in Appendix [X| 



II. MODEL DEFINITIONS 

The class of models is most simply introduced in terms 
of the growth of a one dimensional interface, in which 
both adsorption and desorption processes take place. In 
the present models the key feature is that desorption may 
take place only at the edge of a plateau, i.e. at a site 
which has at least one neighbor at a lower height. We 
study two particular models in this class, (a) a restricted 
solid on solid (RSOS) version where the height differences 
between neighboring sites are restricted to take values 
1, 0, —1 and (b) an unrestricted model where there is no 
such restriction. Both models are defined on a Id lat- 
tice of N sites i — 1 . . . N and associated with each site 
is an integer height variable hi which may take values 
Q, 1 ... oo. Periodic boundary conditions are imposed so 
that h Tv-j-i = hi. 

The dynamics are defined through the following algo- 
rithm: at each update choose a site i at random and 
carry out one of the following processes: Adsorption of 
an adatom 



hi 



hi + 1 with probability q 



(1) 



or desorption of an adatom from the edge of a step 



hi — > min(/ii, /li+i) with probability (1 — q)/2 (2) 
hi — > min(/ii, hi—i) with probability (1 — ?)/2 (3) 

The update of the chosen site i is conveniently imple- 
mented in a simulation by drawing a random number 
between and 1 from a flat distribution. If the number 
is less than q process (Q) is executed, if the number is 
greater than (1 + q)/2 process (||) is executed, otherwise 
process (||) is executed. 

In the RSOS version, a process is only carried out if it 
respects the constraint 



< 1 



(4) 



For both models, when the growth rate q is low, the 
desorption processes (||)-(||) dominate. If all the heights 
are initially set to the same value this layer will remain 
the bottom layer of the interface. Small islands will grow 
on top of the bottom layer through the process (0) but 
will quickly be eliminated by desorption at the island 
edges. Thus, the interface is effectively anchored to its 
bottom layer and the velocity, defined as the rate of in- 
crease of the minimum height of the interface, is zero in 
the thermodynamic limit. On a finite system, very large 
fluctuations will occasionally occur which allow a new 
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FIG. 1. (a) Typical configuration of the interface, nk is 
the fraction of sites at height k above the minimal height in 
the configuration (here no = ni = |, ri2 = The average 
island size grown on top of level k is Ik- (b) Mapping of the 
configuration of Fig. la to the charged-particle representation 
(see Eq. [y]), along with a site coloring, as described in Sec. [v| 



layer to be completed and the velocity will be positive 
but exponentially small in the system size. 

As q is increased the production of islands on top of the 
bottom layer increases until above q Cl the critical value 
of q, the islands merge and new layers are formed at a 
finite rate giving rise to a finite interface velocity in the 
thermodynamic limit. Thus, above the transition one 
expects the velocity to behave as 



v ~ (q - q c ) y , 



(5) 



where y is the critical exponent describing the growth in 
velocity. 

Another critical exponent is defined by considering no, 
the fraction of sites at the lowest level (see Fig. |l|) . Below 
the transition (q < q c ) this fraction will be large since 
the interface is anchored at this level. As the transition 
is approached more and more islands form on top of the 
bottom layer and the fraction no will decrease until it 
vanishes at the transition. This may be described by 



n ~ (q c - qf 



(G) 



Similarly one may define a family of exponents {xk} 
describing the vanishing of n^, the fraction of sites at 
level fc, 



n h ~ (<? c - q) x " ■ 
The interface width is defined by 

W = [N~ 1 '£(h i -N- 1 '£h j ) 



2ll/2 



(7) 



(8) 



Below the transition the width should be finite indicating 
a smooth interface exploring only a finite height above 
the bottom layer. However, above the transition the in- 
terface should display the behavior generic to moving in- 
terfaces Q , that is roughening where the width diverges 
as 



W ~ N 1 ' 2 atq>q c . 



(9) 



Therefore near to and above the transition we expect 

W -N^iq-qJ* , (10) 

where x is the exponent describing the vanishing of the 
roughness as the transition is approached. 

The above considerations hold for both RSOS and un- 
restricted models and and we will address the question as 
to what extent the two models can be seen as represen- 
tatives of a whole class of models with the same critical 
behavior. In Sec. IV we will further explore the scal- 



ing behavior, adding to the exponents y,Xk,X that we 
have so far introduced. However, we defer this until after 



Sec. Ill where we discus the relation to a directed per- 
colation model through which some of the simple critical 
behavior may be understood. 

For the moment we note that the RSOS version (Q) 
may be viewed as a driven diffusion model of two oppo- 
sitely charged types of particles. The charges 



= h. 



e {- , o , +} 



(ii) 



are bond variables and represent a change of height be- 
tween adjacent interface sites (see Fig. [I]). In this repre- 
sentation, it is convenient to convert the dynamical rules 
(§])-(|3|) into processes occurring at bonds with the follow- 
ing rates 



+ — * +0 with rate g 

+ 0^0+ with rate 1 

— 0^0— with rate g 

— — * — with rate 1 

— * H with rate g 

1-^00 with rate 2 
h — > with rate g 



(12) 

where the growth rate g is related to q of equations @)- 
(|) through 



g = 2q/(l - q) . 



(13) 



The rule that desorption cannot occur from the middle 
of a plateau corresponds to the absence of the process 
0^ h In this version the dynamics may be per- 
formed without reference to the height variables hi (al- 
though these can easily be reconstructed to within the 
height of the bottom layer from the variables Ci). The 
critical behavior is reflected in the correlations between 
the charges. At q < q c , the charges are arranged as 
closely bound H dipoles. At q > q c , the dipoles be- 
come unbound, and the fluctuations in the total charge, 
measured over a distance of order N, diverge with N. 

The charged particle representation also allows the ef- 
fect of desorption from the middle of a plateau to be 
studied. This is done in Appendix A by introducing a 
non-zero rate p for the process — ► — + and solving 
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the steady state exactly in the case p = 1 — g/2. The re- 
sult is that, although different choices of this rate allow 
the interface velocity to be positive, zero or negative, the 
interface is always rough and no smooth phase exists. 



III. RELATION TO DIRECTED PERCOLATION 
AND CONTACT PROCESSES 

Some of the critical behavior described in the previous 
section may be related to that of directed percolation 
(DP). In DP sites of a lattice are either occupied by a 
particle or empty. The dynamical processes are that a 
particle can self-annihilate or produce an offspring at a 
neighboring empty site. If the rate of self-annihilation is 
sufficiently high the system always reaches an absorbing 
state where no particles remain and therefore no further 
particle can be produced. On the other hand when the 
rate of offspring production is high, another steady state 
of the system, where the density of particles is finite, ex- 
ists on the infinite lattice and is termed the active phase. 

In DP the dynamics is usually carried out in paral- 
lel, e.g. pi] ]. In the corresponding models in the math- 
ematical literature, known as contact processes Jl9[] , the 
dynamics are defined in continuous time which can be 
numerically implemented by random sequential dynam- 
ics. 

To see the analogy with the growth model defined in 
Sec. H consider the dynamics of the bottom layer of the 
unrestricted model. Let us introduce variables {si} which 
take value 1 if the height of site i is that of the bottom 
layer and take value zero otherwise. The algorithm stated 
in (|l])-(|j) is exactly equivalent to the following dynamics 
for the {si}. At each update randomly select a site i and 
modify Si with the following probabilities 

if Si = 1 

Si — > with prob. q 
if {si-i,Si,s i+ i} = {0,0, 1} 

Sl -> 1 with prob. (1 - q)/2 (14) 
if {si_i,Si,s i+ i} = {1,0,0} 

Si — > 1 with prob. (1 — q)/2 
if {sj_i,Sj,Sj+i} = {1,0, 1} 

Sj — > 1 with prob. 1 — q 

This dynamics is exactly that of a contact process when 
we take Sj = 1 to indicate the occupation of a site: the 
particles self-annihilate with rate A = q/(l — q) and a 
particle is created at a vacant site with rate 1/2 if one 
neighbor is occupied and rate 1 if both neighbors are 
occupied. This process has been extensively studied by 
series expansions |25| and short time expansions |2(| and 
the transition found to occur at A c ~ 0.3032 correspond- 
ing to q c ~ 0.2327 for the unrestricted growth model. 
Thus, as seen by the bottom layer of the growing in- 
terface, the transition from anchored to moving phase is 



simply a DP transition. The anchored phase corresponds 
to the active DP phase whereas the moving phase corre- 
sponds to the absorbing DP phase. 

The critical behavior of DP may be described as fol- 
lows. Above the transition (A > A c ) an initial seed parti- 
cle will typically produce activity over a region of lateral 
extent £j_ ~ H ^ and duration ^ ~ |e| , before the 
absorbing state is reached. Here e, given by 

e = 9-9c, (15) 

measures how far the system is from criticality and £j_ 
and £|| are interpreted as spatial and temporal scaling 
lengths which diverge at the transition. Below the tran- 
sition (A < A c ) the density n of occupied sites in the 
active phase is n ~ |e|^ and there is a finite probability 
\e\P that an initial seed particle will result in the active 
phase being reached. The lateral extent of such an ac- 
tive cluster grows with time as t z where z = V\\jvy_ is the 
dynamic exponent. The typical size I of regions contain- 
ing no activity diverges as the transition is approached 
as I ~ n" 1 ~ \e\-P, 

These exponents allow one to readily identify y and Xq 
defined in Eqs. @ and (||). Since no corresponds to 
the density of occupied sites in the DP active phase we 
expect in (@) 



no ~ (q c - q) 13 and x = j3 . 



(16) 



In the moving phase the velocity is given by the inverse of 
the typical time for the bottom layer to be covered over. 
We identify this time scale with the lifetime of active 
regions in the DP absorbing phase and we expect in (||) 



~ (9 - IcT 11 and y = u\\ 



(17) 



Thus we see that the two exponents xq and y that de- 
scribe quantities involving only the dynamics of the bot- 
tom layer may be directly identified with known (numer- 
ically) DP exponents. Exponents describing quantities 
involving higher levels of the growth process such as Xk 
given by ([?]) are not so straightforward as we shall see in 
Sees. [TV] and [v| 

In this section the exact mapping of the bottom layer of 
the unrestricted growth model to DP has been described. 
For the RSOS model there is no such exact mapping, nev- 
ertheless we expect the bottom layer to exhibit DP be- 
havior and relations (|l^),(|l7|) to hold, because the phase 
transition in this model should display a robustness with 
respect to the microscopic rules similar to that found in 
DP models. 



IV. SCALING AND NUMERICAL RESULTS 

So far we have seen that as criticality is approached 
from the smooth phase the scaling properties of some 
quantities involving only the bottom level of the inter- 
face may be adequately described using DP exponents. 



4 



However, for more general quantities the scaling is less 
clear cut and indeed only a partial picture emerges. We 
first deal with properties of the first few layers in the 
smooth phase using heuristic arguments and then com- 
paring them to numerical results. The width W, a quan- 
tity involving all layers, is studied next. It is argued 
that a simple scaling picture, involving only the DP cor- 
relation length does not hold. We provide evidence 
to suggest that other length scales may be important as 
criticality is approached from within the moving phase. 
A partial scaling picture which emerges is then summa- 
rized. 



A. Scaling properties of the first few layers 

We now discuss the scaling properties of the first few 
layers k — 1,2,... above the bottom layer. Since, in the 
smooth phase, the scaling properties at the bottom layer 
(k = 0) are completely characterized by the three DP ex- 
ponents xq — /3, i>j_, and ui\, it is natural to assume that 
the next layers obey similar scaling laws with analogous 



exponents, Xk — , v 



(V) 



,(k) 



where, for example, 



Xk is the density exponent defined in Eq. (^). In principle 
all these exponents could be different and independent. 
Our numerical results, however, suggest that the scaling 
exponents v± and are identical on all levels and 
equal to the DP exponents v± and v\ | . This remarkable 
property implies that the growth process at criticality is 
characterized by a single anisotropy exponent z — vu/v±. 
On the other hand we find numerically that the density 
exponents Xk for k — 1,2,... are different and consid- 
erably reduced compared to their DP value xq = f3. It 
appears that these exponents are non-trivial in the sense 
that they are not simply related to DP exponents. 

In order to explain the reduced values of Xk , we present 
a simple heuristic argument [ ji7| . Let us first consider 
the bottom layer k — 0. As explained in the previous 
section, DP is characterized by two length scales: the 



and 



average size of inactive islands Iq ~ n^T r 
the spatial scaling length £j_ ~ \t\~ w± They both 

are related, for a system of size N, by the finite size 



-P 



scaling relation Iq ^ 
satisfying f(s) ~ s /3 ^ u± for 



- ) , where / is a function 
0, implying that on a 



finite system at criticality one has uq ~ l 



]\[-fJ/v±_ 



Now consider the next layer k = 1. One may view the 
sites at heights k > 2 as islands of active sites growing 
on top of the inactive islands of the k = 1 level, whose 
typical size is Iq. Applying the same scaling relations 
and assuming that the system size may be replaced by 



Iq, we find n\ 



I, 



where l\ is the mean 



size of islands of sites with height k > 2. Repeating this 



argument one obtains rik 



Ik 1 



-P/ux 
k-l 



and therefore 



n k ~\e\ x \ x k =mivi_t . 



(18) 




(R) 



100 

/ 



k=0 




: k=1 




- k=2 


(U) 


* k=3 









: ° 00 °°°o„ 




k=0 




k=1 




< k=2 


(R) 


k=3 





10 100 1000 10 100 1000 

t t 
FIG. 2. Dynamic simulations: the densities and the 
integrated densities m,k at criticality as a function of time t 
in thefU) unrestricted and (R) restricted growth model. See 
Eq. @). 



Inserting the numerically known DP exponents p8[ 

(3 = 0.27649(4), v\\ = 1.73383(2), v x = 1.09684(2) 

we obtain the approximations x\ ~ 0.070, X2 — 0.017, 
and X3 ~ 0.004 which are in qualitative agreement with 
the numerical results (see below). However, this simple 
scaling argument is not quantitatively correct for several 
reasons. First we consider the temporal average of the 
island size Iq as a fixed 'finite size' length for a DP pro- 
cess taking place on top of the island; we thus neglect 
the temporal fluctuations of Iq. In addition, unlike ex- 
act scaling relations, our scaling argument is expected to 
fail in higher dimensions since it is derived in the case of 
one dimension where active sites separate inactive island. 
Nevertheless it can be used as a rough approximation as 
well as qualitative explanation for the strongly reduced 
values of x k compared to xq for k > 1. 



B. The first few layers - numerical results 



In order to determine the density exponents Xk and 



to verify that the scaling exponents emu ^ 
indeed identical, we employ three different variants of 
Monte-Carlo (MC) simulations, termed dynamic, static 
and finite size simulations, described as follows. 

Dynamic simulations: First we measure the temporal 
evolution of the densities at criticality in a large sys- 
tem, starting from a flat interface and averaging over a 
large number of runs. For times larger than some tran- 
sient time the densities are expected to decay according 
to 



and 



,0) 



are 
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e e 
FIG. 3. Static simulations: double logarithmic plot of the 

saturation value of the integrated densities rrik vs |e|. See 

Eq. @. 




10 100 1000 10 100 1000 

N N 
FIG. 4. Finite size simulations: Measurement of the inte- 
grated densities rrtfc and the growth velocity v at criticality as 
function of the system size N. The solid (dashed) lines refer 
to the unrestricted (restricted) model. See Eq. (fell). 

n k ~ t- x */^l . (19) 

Static simulations: We also determine Xk directly in off- 
critical (static) simulations, measuring the densities n k 
in a sufficiently large system in the smooth phase and 
averaging over very long times. Using this method we 
can determine the exponents Xk directly through the ex- 
pected behavior 

n k ~ |ef * , (20) 
where e = q — q c < 0. 

Finite size simulations: Finally we analyze the finite size 
scaling of n k in critical systems of size N averaged over 
long times. Here the expected scaling behavior reads 

n k ~ N- X ^ v± . (21) 

Thus the dynamic simulation should yield a numeri- 
cal value for x k /v\\, the static simulation a value for x k 
and the finite size simulation a value for x k /v±. If, on 
inserting the DP values of v\ \ and v± , the three methods 
(|l9|)-(^l[) lead to the same numerical result for the ex- 
ponent x k , to within numerical errors, we may conclude 
that the scaling exponents and are indeed equal 
to the DP exponents. 

We observed that the quality of the numerical results 
can be improved considerably if we measure the inte- 
grated densities 







dynamic method 


static method 


finite size method 


x 


U: 


0.275(5) 


0.273(10) 


0.276(5) 




R: 


0.270(10) 


0.277(10) 


0.265(10) 


XX 


U: 


0.114(5) 


0.110(10) 


0.125(5) 




R: 


0.108(10) 


0.110(10) 


0.118(10) 


X2 


U: 


0.039(15) 


0.035(15) 


0.045(10) 




R: 


0.022(15) 


0.025(20) 


0.033(15) 


%3 


U: 


0.011(10) 


no result 


0.015(10) 




R: 


no result 


no result 


no result 



TABLE I. Numerical estimates for the density exponents 
xo,...,X3 for (U) the unrestricted and (R) the restricted 
version of the growth model obtained by various simulation 
methods. 



k 

TOfc = ^ n ( 22 ) 
3=0 

instead of the densities nk ■ The integrated density is 
the probability of finding the interface at height h < k. 
Since Xk-i > Xk, we have nk~i <C near criticality and 
therefore mk and n& scale asymptotically with the same 
exponents. The difference between the two quantities is 
illustrated for the case of dynamic simulations at q c in 
Fig. H (for the method of determining q c see below). As 
one can see, the graph for the integrated densities m*, 
in the unrestricted model shows almost perfect straight 
lines in a double logarithmic representation while the cor- 
responding curves for m, and n% still increase which 
makes it impossible to determine the corresponding ex- 
ponents with the 1000 time steps illustrated. The same 
observation, although with less numerical accuracy, holds 
for the restricted model. Therefore, instead of nk, we al- 
ways measure the integrated densities mk for which we 
assume the same scaling as in Eqs. (|T9|)-(pl]) . 

The dynamic MC simulations are performed on a large 
system with N = 10 4 sites starting from a flat inter- 
face. Detecting deviations from the power law behavior, 
we find the critical points q c = 0.23267(3) for the un- 
restricted and q c = 0.1889(1) for the restricted model 
(corresponding to g = 0.4658 in Eq. (|l3|)). The time 
dependence of mk at criticality is averaged over typically 
10 5 independent runs. The results are shown in Fig. |^. 
From the slopes we estimate the critical exponents Xk/v\\- 
Using of DP we obtain the exponents Xk which are 
summarized in the left column of Table I. 

We note a numerical puzzle we have so far failed 
to explain. The critical value of q for the restricted 
model appears to be given to high numerical accuracy 
by <7u/(l — <7u) where q u is the critical value of the un- 
restricted model. Whether this is sheer coincidence or 
whether it is a manifestation of some duality between 
the two models is an open question. 

Static simulations are then carried out in the smooth 
phase q < q c . Varying |e| — q c — q from 10~ 3 to 10 _1 
we first equilibrate the interface on a large lattice with 
N = 10 4 sites over a time interval up to 10 5 time steps. 
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FIG. 5. Finite time static simulations: Data collapse for 
the scaled integrated densities mt as a function of time mea- 
sured in the smooth phase for various values of e in the case 
of (U) the unrestricted and (R) the restricted growth model. 



Then the stationary densities rrik are averaged over a 
time interval of the same size. The results are again 
averaged over 100 independent runs. Using this method 
we can estimate the critical exponents Xk directly from 
the slopes of the lines in Fig. (see middle column in 
Table I). 

Finally finite-size simulations at criticality are car- 
ried out for various lattice sizes N — 8, 16, 32, . . . , 1024. 
Starting from a fiat interface we averaged the integrated 
densities rrik over time intervals proportional to N z , rang- 
ing from 5000 time steps for N = 8 up to 3 • 10 7 time 
steps for N = 1024. Since finite sized systems at criti- 
cality have a small but finite growth rate, the densities 
mk have to be measured with respect to the actual low- 
est level of the interface, i.e. in a co-moving frame. The 
slopes of the curves in Fig. |^ give an estimate of x^/vx 
(and therewith Xk , see right column of Table I) . In addi- 
tion, the finite growth rate measured in this type of sim- 
ulation is expected to scale as d ~ N~y/ U ± which allows 
an estimate of the exponent y in Eq. (0). Our results 
are y = 1.71(5) for the unrestricted and y = 1.76(10) for 
the restricted model (see Fig. ||) which is in agreement 
with our scaling prediction y = i>\\ ~ 1.734 in Eq. (|l7|). 

As a final check of scaling we can also obtain a finite 
time collapse using the short time data from static sim- 
ulations. In the smooth phase the expected scaling is 
mk{i) ~ |e| Xfc /fc(t|e| I/ ll), where fk(s) — > const for s — > 00 
and jffc(s) ~ s~ Xk / u w for s — > 0. Fig. || shows the scal- 
ing functions /fe(s) for various values of e ranging from 
— 10 -4 to —0.07 measured up to 10 4 time steps. The long 
time stationary values of the densities in static sim- 
ulations ((2^) correspond to the saturation levels of the 
different curves. 

Since the three different methods lead to the same re- 
sults for Xk within numerical errors, we conclude that the 



scaling exponents and v)?' are indeed identical on all 
levels and equal to the DP exponents v±_ and v\ \ . We ob- 
tain the DP exponent xq = /?, as expected at the bottom 
layer, while x\ ~ 0.12, x 2 ^ 0.04, x 3 ~ 0.015, . . .. These 
latter exponents, although showing the same trend, are 
clearly different numerically from those obtained using 
the heuristic picture of subsection IVA. It is not clear 
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FIG. 6. Width at criticality as a function of time in the 
unrestricted model. The graph is for TV = 2048 and shows 
the growth of the width starting from a fiat interface averaged 
over 2000 runs. The straight line is a best fit through the long 
time points (2 6 — 2 13 MCS) and has a slope of 0.24. 



whether {xk} for k > are related to the DP exponents 
or whether they are independent exponents. The fact 
that we obtain the same exponents in the restricted and 
the unrestricted model suggests that both variants - at 
least with respect to the first few layers - belong to the 
same universality class. 



C. Scaling of the Interface Width 

In this subsection we investigate the scaling properties 
of the interface width, defined by (j^) , at criticality and 
as the rough phase is entered. The latter measurements 
lead us to conclude that a complicated scaling behavior 
prevails. 

First recall that in the smooth phase the interface ex- 
plores only a finite number of levels and the width is 
finite. In the rough phase the width is expected to di- 
verge according to W ~ N 1 / 2 , the behavior generic to 
moving interfaces in one dimension || . 

A naive scaling picture would suggest that the inter- 
face width may be written as W ~ \ e \~ n f(N/£±) ~ 
\e\~ v f(N\e\" ± ), where rj is some critical exponent. We re- 
fer to this picture as naive because implicit is the simple 
scaling assumption of a single length scale in the prob- 
lem, a point we shall question in Sec. IV D. Within this 
naive picture one chooses a scaling function / so as to 
obtain the expected v^V behavior in the rough phase. 
This leads to the following small e, large N asymptotic 
limits: 



W ~ |e|-'> for e < 

W - Ni/ u -l for e = 

W - /V 1 /2 e ^/(2- I? ) for e > o 



(23) 
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FIG. 7. Width at criticality as a function of time in the 
restricted model. The graph is as for Fig. ^ and the system 
size is N — 2048. The best fit through the longer time points 
(2 6 — 2 13 MCS) has a slope of 0.43. 



We now proceed to examine the actual behavior of the 
width at criticality. In Figs. width against time 

is plotted for simulations run at criticality. For both 
unrestricted and restricted models a long time scaling 
behavior emerges before the width saturates due to the 
finite length of the system. For times shorter than the 
saturation time, the width at criticality, W Cl grows as 



W c ~ (lnt)" 



(24) 



Now since the saturation time is expected to scale as N z , 
where z is the dynamic exponent, one deduces that the 
saturation scaling of the width is 



Wr 



(In TV) 7 , 



(25) 



where 7 is given by 7 ~ 0.24 w 1/4 for the unrestricted 
model and 7 ~ 0.43 for the restricted model. These nu- 
merical results suggest that the critical exponent 77 in 
Eq. (23) is in fact equal to zero. For the restricted model 
the scaling of the width at criticality is similar to that 
of an unrestricted polynuclear growth model || where 
7 = 1/2. The relationship between the prese nt models 
and PNG models will be discussed in Sec. VII] . However, 



for the unrestricted model the value of 7 is clearly dis- 
tinct from 1/2 which shows that the critical width could 
display non-universal behavior. On the other hand 7 
could be restricted to a finite number of possible values. 
It should also be noted that in prelimin ary simulations 
the value of 7 was erroneously assumed Q . 

We now examine how the saturation width diverges as 
the growth rate is increased and the interface enters the 
moving phase. To do this it is convenient to subtract 
out the critical width and measure AW(e) = W(e) — 
W c - Since W c is negligibly small as compared with W(e) 
for e > one expects W and AW to have the same 
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FIG. 8. Double logarithmic plot of width in the moving 
phase as a function of e. The system size was N — 512 and 
each point is an average over 2000 simulations. The simula- 
tions were started from a fiat interface and allowed to saturate 
over 2 13 MCS. The width was then averaged over another 2 13 
MCS. The straight line corresponds to \ = 0.95. Simula- 
tions of the restricted model yielded a similar behavior and 
estimated exponent value \ = 0.92 



asymptotic behavior. In Fig. 
behavior is 

AW(e) ~ 



it is seen that the scaling 



(26) 



where \ = 0.95(5) for the unrestricted model. 

This result for \ disagrees with the naive scaling pic- 
ture ( p3[ ) which suggests that the width should diverge 
with the exponent vj_/2 = 0.55 (given that 77 was found 
to be equal to zero). We are left to conclude that a more 
complex behavior than simple scaling takes place. 



D. Length Scales in the Moving Phase 

In the previous subsection we saw that a simple scal- 
ing argument involving only the DP scaling length does 
not describe correctly the numerical results. In this sub- 
section the aim is to identify possible additional length 
scales and indicate that complex critical behavior may 
be present. Therefore the subsection is by nature specu- 
lative. 

We investigate scaling lengths of the first few layers 
in the moving phase. We consider the correlation func- 
tions (mfc(i)TOfc(i-|-r)) of the integrated densities ( p^ ) 
between sites i and i + r. 

Let us first consider the bottom layer density mo = no- 
In the smooth phase it has already been noted the dy- 
namics of the bottom layer (of the unrestricted model) 
is exactly that of DP in the active phase. However, in 
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FIG. 9. Correlation functions {rrik(i)mk(i + r)) / {rrik} in 
the moving phase in the unrestricted model. The system pa- 
rameters are N = 2048, e = 0.02 and the results are an aver- 
age over 1000 simulations each run for 2 16 MCS to equilibrate. 



the moving phase there is a subtle difference between the 
present model and DP in the absorbing phase. This dif- 
ference stems from the fact that in the present model 
there is no absorbing state. Instead, whenever a layer 
is completed the next highest layer becomes the bottom 
layer. We effectively move in a frame co-moving with the 
lowest uncompleted layer and relabel the layers appro- 
priately. 

In order to show that this produces a non-trivial ef- 
fect we plot {mk{i)m,k{i + r)) / {rrik) in Fig. ^|. First recall 
that in the DP absorbing phase near criticality one ex- 
pects to see a power law decay (no(i)no(i + r)) / (no) ~ 
(mo(i)mo(i + r)) / (mo) ~ r ~P/ u ± f the correlation func- 
tion up to a scaling length which diverges as e - " 1 . At 
lengths longer than the scaling length the correlation 
function should decay exponentially with r. In Fig. ^ 
for k — we see an initial power law decay with power 
given by —0.27, then a crossover at r ~ £j_ to a new 
power law decay with power ~ —0.76 (rather than to 
an exponential decay as would be the case with usual 
DP scaling). We checked for different system sizes and 
values of e that this qualitative behavior (crossover to a 
new power law at long distance) was reproduced. Also 
on Fig. H correlation functions for the integrated densities 
at higher levels are plotted and again one sees crossovers 
between two power laws. The two powers appear distinct 
for each level k. For k = 1 the two powers are 0.15,0.56 
and and for k = 2 the two powers are 0.09,0.36. It ap- 
pears that the crossovers occur at a length dependent on 
k although it is difficult to quantify this. If the crossover 
lengths were k dependent then it would imply different 
length scales existing on different levels 

In order to give a heuristic explanation for the above 
behavior let us consider the density at the bottom layer 



too = riQ. A picture that could explain the crossover phe- 
nomena in Fig. H is that the sites at the bottom level are 
grouped into clusters. Each cluster displays the scaling 
behavior of an active DP cluster and therefore is of typi- 
cal size e~ Vx and within a cluster the correlation function 
decays as ~ r -P/"± , Thus the correlations within clus- 
ters generate the first power law decay. However, one 
also has correlations between clusters. Thus within the 
cluster picture the second power law measures the decay 
of correlations between clusters. 

A similar cluster picture could hold for the first and 
second layers k = 1,2. Within the picture, the sites at 
the first level, for example, are distributed in clusters. 
The first power law in the k = 1 curve of Fig. ^| measures 
the correlation within a cluster and the second measures 
correlations between clusters. The fact that the crossover 
appears to occur at a different distance r than for the 
bottom layer indicates that the clusters at the first level 
are larger than those at the bottom level. Therefore there 
is more than one length scale in the problem. 

This picture is far from being complete or verified un- 
ambiguously and many questions remain open. For ex- 
ample it is not clear whether the second power laws in 
Fig. ^| continue indefinitely or are cut off at some larger 
length. The numerical value of x is also not explained. 



E. Summary of the Scaling Picture 

For the sake of clarity we summarize the scaling results 
of this section and the partial picture of scaling. First we 
have the smooth phase with exponents Xk associated with 
the density at each level, xq is given by the DP expo- 
nent (5 whereas Xk for k > 1 appear to be non-trivial, in 
the sense that they are not simply relate d to D P expo- 
nents. The simple approximation of Sec. IV A gives the 
qualitative trend but is ruled out quantitatively by the 
numerics. As the roughening transition is approached 
the DP scaling length £j_ and scaling time £|| hold at all 
levels. This implies that approaching the transition the 
dynamic exponent is z = v\\jv^. 

At the transition the interface has logarithmically di- 
verging width of the form (|25|). However the value of 7 
appears to depend on which version of the model one sim- 
ulates. This could either point to 7 being non-universal 
or that it can take one of a finite number of values. 

Above the transition the velocity grows with the DP 
exponent v\\. This reflects presence of the DP scaling 
length and time at the bottom level. However by mea- 
suring the growth of the interface width we have ruled 
out a simple scaling picture involving only the DP scaling 
length. We have provided evidence to suggest that there 
may be longer scaling lengths which characterize the size 
of clusters of sites at higher levels. This picture remains 
to be fully investigated. Likewise the question of scal- 
ing times for different levels in the rough phase has not 
been fully addressed. (In the rough phase the dynamic 
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exponent should take the KPZ value z = 3/2). 

It is interesting to note that in a very recent renormal- 
izati on g roup analysis of the field theory introduced in 



o.oyi 



Sec. VII, multicritical behavior was identified |32)| . This 
could be consistent with the complicated scaling behavior 
we have observed. 



V. SPONTANEOUS SYMMETRY BREAKING 

The growth models defined in Sec. |l|may be viewed as 
examples of spontaneous symmetry breaking (SSB) in a 
Id system with periodic boundary conditions and a non- 
conserved order parameter. As will be shown, the models 
display a robust local mechanism for eliminating islands 
of the minority phase generated by fluctuations. 

The symmetry of the growth models, apart from spa- 
tial translation and reflection invariance, is a (discrete) 
translational invariance in the heights (Zoo). In the 
smooth phase this symmetry is broken since the system 
spontaneously selects one of the heights as a reference 
level which then serves as bottom layer for local fluctua- 
tions of the smooth interface. 



A. Order parameters. 

In order to quantify the symmetry breaking, we define 
a magnetization-like order parameter (valid for both the 
restricted and the unrestricted models): 



Mi 



1 N 

-Y( 

N ^ K 



(27) 



This order parameter is clearly not conserved by the dy- 
namical rules of the models. In the smooth phase q < q c , 
it has a nonvanishing expectation value (Mi) ^ in the 
thermodynamic limit where (. . .) denotes thermal aver- 
age. On the other hand, in the rough phase the interface 
explores many height levels therefore the contribution to 
the magnetization at different sites are uncorrelated over 
long distances and (Mi) = 0. Near the phase transition 
we expect the order parameter to vanish as 



(Mi 



(28) 



where 9\ is the associated critical exponent. 

The order parameter and the SSB mechanism are il- 
lustrated in Fig. |Io|. Here the heights are represented 
by alternating black and white coloring (c.f. Fig. |l|) and 
therefore the average greyscale in the figure is related to 
the magnitude of M\, Also it is shown how a large is- 
land of one phase shrinks when introduced into a system 
dominated by the other phase, thus ensuring the stability 
of the smooth phase. This behavior is typical of islands 
of all sizes, except for those extending over the whole 
system. 
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FIG. 10. MC simulation of the restricted model for a sys- 
tems of size 600 at different values of q. Each configuration is 
a row of pixels, with sites at even (odd) heights represented 
by black (white) pixels, visualizing the order parameter M\. 
Configurations at intervals of 7 moves per site are shown up to 
2100 sweeps. Initially, a large island of size 400 is introduced. 
For q < q c , the island shrinks and disappears, illustrating the 
mechanism that insures long-range order in the smooth phase. 
The typical time it takes for the island to disappear depends 
on q, it increases and finally diverges when q — > q c . Simi- 
larly the magnitude of the order parameter Mi (visualized by 
the greyscale contrast between the two phases) decreases. At 
criticality the order parameter is zero so that the island is not 
visible any more. 
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0.64(3) 
0.66(6) 


0.40(2) 
0.37(4) 


0.24(2) 
0.21(4) 


0.15(1) 
0.14(3) 


0.11(1) 
0.10(2) 



TABLE II. Numerical estimates for the order parameter 
exponents 8^ obtained from dynamic MC simulations. 



Since the symmetry of the model in the heights is , 
we can define a family of order parameters which gener- 
alize Mi: 



M„ 



1 /2nihi 

3 = 1 



n+1 



(29) 



These order parameters have the same qualitative be- 
havior as Mi and can be understood as discrete Fourier 
transforms of the height probability distribution. It turns 
out that each order parameter is characterized by a dif- 
ferent critical exponent: 



(M n ) 



(30) 



As in the case of the density exponents Xk, we deter- 
mined the exponents 9 n numerically by static, finite-size 
and dynamic MC simulations. The most precise data are 
obtained from dynamic simulations (see Fig. [ll|). The 
numerical estimates for 9 k are summarized in Table II. It 
seems that these exponents are non-trivial in the sense 
that they cannot be related in a simple manner by scaling 
relations to the known DP exponents. 
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FIG. 11. Measurement of the order parameters (M n ) in 
dynamic MC simulations for (U) the unrestricted and (R) the 
restricted variant of the growth model. 

B. External ordering fields. 

For each order parameter M n one can define a conju- 
gate ordering field that favors states where the order pa- 
rameter is positive. This ordering field can be introduced 
by a periodic modulation of the growth rate, i.e. we re- 
place the uniform growth rate q by a height-dependent 
growth rate 

q^q(hi) = q-Xcos(^p^), (31) 
V n + 1 / 

where A is the magnitude of the ordering field. For ex- 
ample, for n = 1 and A > the growth on the bottom 
layer and other even layers is penalised whereas growth 
on odd layers is favored. At criticality, the response to 
this external field is expected to obey a power-law be- 
havior 

M n ~ A K " . (32) 

We measured the exponents K n in static MC simulation 
at criticality (see Fig. |l2|). Varying A over two decades 
from 10~ 3 to 1CP 1 wc obtain the estimates k± — 0.60(4), 
ki = 0.42(3), m = 0.26(3), K X = 0.17(3), and Ki = 
0.12(3). Comparision with the results in Table II suggest 
that K n n . 



C. Spontaneous breaking of continuous symmetries. 

So far we have shown that the growth models discussed 
in this paper exhibit spontaneous breaking of a discrete 
symmetry in one dimension. It is therefore of interest to 
ask whether a continuous symmetry can also be broken 
in one-dimensional nonequilibrium models. Continuous 
symmetries, where the order parameter can take a con- 
tinuum of values, seem to be harder to break than dis- 
crete symmetries. Consider for example the equilibrium 
case: discrete symmetries can be broken above one di- 
mension; continuous symmetries, however, can be broken 
only above two dimensions, although weaker vortex or- 
dering transitions are possible in 2d (Kosterlitz-Thouless 
transitions ps[]). 




X 



FIG. 12. Response of the order parameter Mi to an ex- 
ternal ordering field at the critical point in the unrestricted 
model with 2000 sites. Measuring the slope of the line in the 
double-logarithmic plot over two decades gives estimates for 
the response exponents K n (see text). The saturation for very 
low values of A seems to be a finite-size effect. 

Here we present a version of the growth model with a 
continuous symmetry which is broken spontaneously in 
the smooth phase. The only difference between this ver- 
sion and the unrestricted model described above is that 
the height increment at a growth move is a continuous 
rather than a discrete number. 

The model is defined on a Id lattice with periodic 
boundary conditions and continuous (real) height vari- 
ables hi at sites i = 1 . . . N. The interface evolves by 
choosing a site i at random and carrying out one of the 
processes 

hi — > hi + £ with probability q 
hi — > vain(hi,hi+i) with probability (1 — q)/2 (33) 
hi — ► min(/ii, hi-%) with probability (1 — q)/2 

where £ is a positive real random number selected from a 
flat distribution between and 1. The symmetry of this 
model (apart from spatial translations and reflections) 
is continuous translational invariancc in the heights, i.e. 
overall shifts of the interface heights by any amount. The 
symmetry breaking corresponds to the fact that in the 
smooth phase the interface selects a given height which 
is a real number, and remains pinned to that height for 
a time that grows exponentially with the system size. 

Starting from a flat interface at height zero, the dy- 
namics taking place at the lowest level in the continuous 
model are exactly the same as in the unrestricted version 
of the discrete model, the simple reason being that in 
both cases each height level is decoupled from the higher 
ones. This means that both models have the same criti- 
cal point q c = 0.23267(3). Moreover, they have the same 
occupation of the lowest level which is characterized by 
n = (q c — q) x ° where x — (3 ~ 0.277 is the DP density 
exponent. 
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FIG. 13. Spontaneous breaking of a continuous symmetry: 
The graphs show the occupation of the bottom layer no and 
the growth velocity v in a system of size N at criticality. 
From the slopes we obtain the exponents xo/v± = 0.245(10) 
and y/u± = 2.59(7). 



Above the critical point, the interface is rough and 
has a finite growth velocity. Simulations show that the 
roughness exponent characterizing the interface at this 
phase is consistent with the KPZ exponents 0, within 
numerical accuracy. However, an interesting difference 
from the discrete model occurs in the growth velocity v. 
As in the previous case, it is related to the inverse of the 
island lifetime. However in the present case, the step size 
by which the interface grows is not 1, but rather a real 
random number between and 1. Each completed layer 
corresponds to the growth of the interface by the smallest 
surviving step. Now, from our results on the discrete 
models we expect that at criticality the width behaves as 
some power of In TV. Since the number of steps in a finite 
system is of the order of A', we expect the smallest of 
them to be of order l/N. At the critical point, therefore, 
we expect the following finite size scaling for the velocity 



(34) 



which corresponds to the scaling 

v~{q-q c ) v , y = v\\ +is±_ ~ 2.83. (35) 

This has to be compared with y — v\\ ~ 1.73 in the 
discrete model. These exponents are in good agreement 
with the values measured in MC simulations, as shown 
in Fig. [ll]. It would be interesting to study order param- 
eters M u which are continuous generalizations of M n in 
Eq. (§): 



1 N 

— 22 exp(27riw hj) 



As before, we expect a power law behavior (M^ 
where 0(u>) depends continuously on u>. 



(36) 



|0(uO 



VI. MEAN FIELD APPROXIMATION 

In this section we derive the mean field equations cor- 
responding to the model (0)-(||) and study the resulting 



steady state distribution. We consider the unrestricted 
SOS model since the equations for this case are some- 
what simpler. However, both the restricted and the un- 
restricted models are expected to exhibit the same qual- 
itative dynamical behavior. 

To derive the mean field equations we introduce at each 
site i a set of variables ipk(i), k = 0, ...,oo. Here ipk(i) is 
equal to 1 if the interface at site i is at height k and it is 
equal to otherwise. Let (V'fc(i)) be the average of ipk{i) 
over all realizations of the dynamical equations starting 
with the same initial configuration. Let us first consider 
the occupation of the zeroth level. The adsorption and 
desorption processes defined in Eqs. (0)-(||) result in the 
following equation for (ipo(i))'- 



A 



dt 



Q#oW) 
1 



(37) 



+ ^(l-?)<^)(t)(l-'0o(* + l)), 

where A is a time constant, which, for simplicity, may 
be taken as 1 — q. Within the mean field approximation, 
one replaces the correlation function (i(>o(i)ipo(j)) by the 
product (ipo(i))(ipo(j)) . For the ring geometry considered 
in this work (ip (i)) is independent of the site index i. 
Denoting (ipo(i)) by po we obtain the following dynamical 
equation for po 



dpo 
dt 



-rrr = ~QPo +M 1 _ Po) 



(38) 



where q = q/(l — q). Similar considerations yield the 
following set of equations for pk = ( - 0fc(*)) for k > 1: 



dpk 
dt 



fc-i 



qp k -! - qp k + pk 1 - ^2Pj) ~Pk^2 p i 



(39) 



It is convenient to rewrite Eqs. (38),( |39|) in terms of the 
integrated density variables (c.f. Eq. (122)) 



fc 

<t>k = Y.Pi- 

j=o 

Substituting (|4(]) and using pk = <fik 
(p8|) , (pSj) into the form 



(40) 

b k -i brings 



d(j>o 
dt 

d(j)k 
dt 



£00 



(41) 



e(f> k - 4> 2 k + (1 - e)0 fc _i (fc > 1) 



where e = 1 — q. These equations have a stationary 
solution corresponding to a smooth interface for e > 0. 
The roughening transition takes place at e = 0. The 
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FIG. 14. Mean field approximation of the integratf 
height density (j>k at level k. 
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FIG. 15. Mean field approximation of the magnetization 
Mi as a function of ln[— m(e)]. 



stationary solution for e > may be calculated by tl 
following recursion relation: 



l r 



+ v/e 2 + 4(1 - e)0 fc _! 



(4: 



with 4>o = e. To leading order in e Eq. ( |42| ) takes tl 
form 



<Pk = \/ 4>k-i ■ 



(43) 



The steady state distribution corresponding to this re- 
cursion relation is 



4>k 



- .(1/2)* 



(44) 



Therefore the mean field values of the exponents Xk de- 
fined in Eq. (Q) are given by 



1 

2*' 



(45) 



The integrated height density <pk is a monotonically 
increasing function of k, varying from e for k = to 1 for 
k — > oo (see Fig. |l4|). It exhibits a rapid increase near 
k ~ fc m which is determined by the following equation: 



&& m -i + 4>k m +\ - %4>k m = o . 



(46) 



The index k m corresponds to the interface height at 
which the density pk — 4>k — 4>k—i is maximal. Using 
Eq. (|43|) one finds that to leading order in e 



k m ~ 7~Z M - ln ( e ) 

In 2 



(47) 



It is easy to see that the width of the (f> distribution re- 
mains finite even in the limit e — * 0. The interval Ak 
corresponding to a change of <p from some value ^min 1° 
another, say, 4> ma x is given to leading order in e by 



Afc = -!— ln[ln(^ mi „)/ln((/) TOax )]. 



(48) 



This interval is independent of e, and thus the width of 
the 4>k distribution is finite. This feature, together with 
Eq. ([l?]) yield the following behavior for the average 
height : 



(M~#> ~ln[-hi(e)] 



(49) 



We now turn to the magnetization-like order parameter 
M i = EjC-lVft (see Eq. @). For small e > Eq. @) 
yields the following expression for M% 



Mi = e - 



y^ e (i/2) f 
fc=i 



(1 



c(l/2)' 



)(-l) A 



(50) 



In Fig. [l5| we plot the magnetization as a function of 
ln(— ln(e)). It is readily seen that Mi does not decay to 
zero for small e. Rather it oscillates with an amplitude 
which remains finite at small e. This may easily be under- 
stood in the following way: the main contribution to the 



sum (50) comes from a finite number of layers centered 
around k = k m . As e decreases, k m increases, and the 
magnetization order parameter probes different layers. 
Since the contributions of the layers have alternate signs, 
the magnetization becomes oscillatory. This feature is 
an artifact of the mean field approximation. In this ap- 
proximation fluctuations of the position of the interface 
are neglected. When they are taken into account they 
are expected to broaden the interface, resulting in a di- 
verging width at the roughening transition. This should 
lead to a vanishing magnetization at the transition. In a 
similar fashion, the mean field theory of KPZ interfaces 
does not describe the roughening behavior. 



VII. FIELD THEORY 

In order to understand the universal properties ob- 
served in the growth models at criticality, it is useful 
to study the corresponding field theory. As will be ex- 
plained below, such a field theory describes a hierarchy 
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of unidirectionally coupled DP processes. Thus it should 
play a role in even more general contexts, namely when- 
ever DP-like processes are coupled unidirectionally with- 
out feedback. 

Unlike the KPZ equation Q, the field theory we con- 
sider involves separate fields for each height level in order 
to incorporate the rule that atoms cannot desorb from the 
middle of an island. Let us first try to guess the Langevin 
equations by adding appropriate diffusion and noise con- 
tributions such that (a) the equation for the lowest level 
reproduces the ordinary Langevin equation for directed 
percolation |}0| , (b) at a given time and position the sum 
over all height densities J2 k P k 1S equal to one, and (c) the 
hierarchy of equations is translational invariant in space 
and time as well as in the heights. The simplest set of 
Langevin equations that meets these requirements reads 
(suppressing the arguments x, t) 

d t po = -qpa + Po(l - Po) + V 2 p + Vo (51) 

k k-1 

d t pk = qpk-i - qpk + pk ( 1 - Pi) - Pk^Pi 



i=0 



i=0 



+V 2 Pk - V 2 /7fc-i + r) k - V k-i (k > 1) , (52) 

where n k (x,t) are field-dependent Gaussian noise fields 
with two point correlations to be specified below. Notice 
that Eq. ( pl| ) is just the ordinary Langevin equation for 
DP |[o| . One can also verify that the sum over all density 
fields X)fe=o P k IS a constant of motion. 

As in the previous section, these equations can be sim- 
plified by introducing integrated density fields 4> k {x, t) = 

J2j=oPj( x ^)' resulting in 
d t 4>k = acj> k - 4>l + #fc-i + V 2 (0 fc - 4>k-i) + Vk , (53) 

where a = 1 — q and = 0. Interestingly, the in- 
troduction of integrated densities also led to a consider- 
able improvement of the numerical results in Sect. IV, 
suggesting that these quantities are more natural in the 
context of the present problem than the height densities 
Pk(x,t). 

Although the Langevin equations (p3| ) follow quite nat- 
urally from the principles (a)-(c) stated above, it can 
be dangerous to conjecture the correlations of the noise 
fields r]k(x,t). A systematic derivation of the Langevin 
equations and the noise correlations based on a bosonic 
operator formalism l3^] will be presented in Ref. [§2| 



Dropping irrelevant operators (like V 



>k-i 



in Eq. 



and introducing independent coefficients for all terms it 
is shown that the Langevin equations are given by 

d t (j> k = a k cj) k - \ k <j) 2 k + q k (f> k -i + D k V 2 cf> k + r\ k (54) 
with noise correlations 

(r)k(x,t)r)i(x',t')) = 2T ki i (p k (x,t) 6{x - x') 6(t - t 1 ) , 

(55) 



where k < I. Notice that there are noise correlations 
between different height levels k,l which are generated 
in a one-loop approximation by non-trivial mixed cubic 
vertices p2|. 

In Sect. |IV| we observed numerically that the scaling 
exponents V\_ and v\\ are identical on all height levels. 
This observation can be verified easily within an 'im- 
proved mean field approximation' as follows. Consider 
the scaling transformation 



Ax, t->A z t, cj> k A- Xh <f> k : 



(56) 



where z is the dynamical exponent and Xfc are the scaling 
exponents of the fields 4> k . Under rescaling, Eqs. ( p4[ ) and 
(p5| ) turn into 

d t <f> k = a k A z <f )k -\ k A z -xx<f>t (57) 
+%A*+**-**-^ fe _ 1 + D k A z - 2 V 2 k + v ' k , 

( Vk (x,t)r}l(x',t'))=2A z+ ™- d r k , l( l> k (x,t) (58) 
x S(x - x') 5(t - 1') , (k<l) 

where d is the spatial dimension. Thus an infinitesimal 
rescaling by A = 1 + m would result in a change of the 
coefficients 

d m a k = a k z 

d m Afc = A fc (z - Xk) 

d m D k =D k (z- 2) 

dm r fe ,i = T k ,i (z + xi- d) (59) 
d m q~k = q~k(z + Xk - Xk-i) 

At the critical dimension d = d c , we expect the coef- 
ficients to be invariant under rescaling. As usual, the 
DP equation at the lowest level k = yields the solu- 
tion z = 2, xo = 2, and d c = 4. The linear term is 
relevant so that the parameter ao has to be tuned to 
zero (this is the mean- field critical point of DP). Also 
at higher levels k > the linear term is the most rele- 
vant contribution wherefore a k = is the multicritical 
point in mean-field. Requiring that the nonlinearity and 
the coupling to the previous level are equally relevant 
(z + Xk — Xk-i = z — Xk), we are led to the solution 
Xk = 2 1_fe . Identifying z = v\\/v±_ and Xk = x k /v±_ we 
obtain, in agreement with Sec. VI, the mean field expo- 
nents 



i/jf F = 1 



v MF = - 



MF 



= i ■ (60) 



The field theory (p4[)-(p5|) may be interpreted as a hier- 
archy of DP processes which are unidirectionally coupled 
by the term q k <p k ^i . Identifying the fields </>o, <j>i, 4>2, ■ ■ ■ 
with densities of particles A, B,C, . . ., this field theory 
corresponds to the reaction-diffusion process 

A <-> 2 A , A—>B 
B^2B, B C 
C^2C, C^D 
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MC simulations indicate that this reaction-diffusion pro- 
cess belongs indeed to the same universality class as the 
roughening transition discussed in this paper. We there- 
fore expect that this field theory describes not only the 
present growth models but any system in which DP pro- 
cesses are coupled in one direction. 



VIII. RELATION TO POLYNUCLEAR GROWTH 
MODELS 

In the previous sections we have seen that the dynam- 
ics of the lowest layer undergoes a DP transition that 
corresponds to a transition from zero to finite velocity of 
the interface. For a class of models termed Polynuclear 
Growth (PNG) @,@,|,|l5],[l6|, which employ parallel dy- 
namics, a similar scenario pertains to the growth at the 
highest level. In these models the use of parallel dynam- 
ics implies that the maximum velocity is 1 i.e. the sites 
at the highest level h = T grow at every time step T. 
The sites at the highest level may be considered as ac- 
tive sites of a DP process and below the transition there 
is a non-zero density of such sites. Above the transition 
there are no sites at the highest level and the velocity is 
less than 1. This transition is lost, however, when the 
dynamics are performed random-sequentially since then 
there is no maximum velocity. This contrasts with the 
model (0)~(^) where a phase transition is found whether 
the dynamics be implemented random-sequentially or in 
parallel. 

In this section we shed some light on the connection 
between the present models and models of the PNG class. 
In order to do this we first generalize the dynamics of the 
unrestricted model (0)-(||) to encompass both random 
sequential and parallel dynamics. In a time step all sites 
are updated according to the following rule: 

hi{T + l) 

= hi(T) with prob. 1-A (61) 
= ft,-(T) + l with prob. qA 

= min [hi(T),h l+1 (T)} with prob. (1 - q)A/2 

= imn[h i - 1 (T),h i (T)] with prob. (1 - q)A/2 

As A varies from to 1, the rule (|6l]) interpolates be- 
tween random sequential dynamics and parallel dynam- 
ics: for A <C 1/N the height of at most one site is modi- 
fied at any time step and the dynamics becomes random 
sequential; for A = 1 all sites are modified at each time 
step and the dynamics is parallel. 

We now make a transformation suggested to us by J. 
Krug (private communication). If one defines 

hi(T)=T- gi (T), (62) 

the variables gi{T) undergo the following dynamics: 

9i(T + l) 



= g l (T) + 1 with prob. 1 - A (63) 
= gi(T) with prob. qA 

= max[ gi (T),g i+1 (T)] + 1 with prob. (1 - q)A/2 

= max [cfe_i(T),pi(T)] + 1 with prob. (1 - q)A/2 . 

The rules (|63] ) yield a growth model where the maxi- 
mum velocity is 1 and the DP transition appears at the 
maximal height level, as occurs with PNG models. How- 
ever, it is important to note that the dynamics j6^ ) is 
always parallel in nature since the heights of many sites 
are modified at each time step. Thus both random se- 
quential and parallel versions of ( |6l| ) are mapped onto a 
parallel rule (]63]) . 



We now compare (63) with a specific PNG model stud- 
ied in ||. In that model the heights of a Id interface are 
updated in parallel in two sub-steps. First all up (down) 
steps of the interface move deterministically to left (right) 
a distance of u lattice spacings. Then all heights are in- 
cremented by 1 with probability p. (For our purposes 
it is convenient to group the two substeps in the reverse 
order to that of || but this does not change the dynam- 
ics). Thus the model is unrestricted. For the case u = 1 
this dynamics may be written as a single parallel update 
where 

9i(T+l) (64) 
= max [pi_i(T), gi(T),g i+1 (T)} + 1 with prob. p 
= max Lg ( -i(T), gi(T),g i+ x(T)] with prob. (1 - p) . 



Clearly the two dynamics (63),(|64j) are similar in char- 
acter but distinct. Furthermore in || it was found 
that at p c ~ 0.527 the width of the interface scales as 
W ~ (In A) 1 / 2 which is distinct from the behavior of the 
unrestricted model but similar to that of the restricted 
model of the present paper. 

In this section we have seen a subtle connection be- 
tween the growth models studied in the present paper 
and models similar in character to PNG models. This 
was done by generalizing the dynamics of the present 
model to a dynamics that interpolates between random 
sequential and parallel. Then, transforming the present 
model via ( |62|) one obtains a model with parallel dynam- 
ics which has similarities with, but is distinct from, the 
PNG models previously studied. However, with the stan- 
dard PNG model ( |64] ) the dynamics is strictly parallel 
and it is not clear if it can be transformed to any ran- 
dom sequential model. Indeed, simply employing the in- 
verse transform of (^32j) yields another model with strictly 
parallel dynamics. 



IX. CONCLUSIONS 

In summary, we have studied in detail a one- 
dimensional stochastic growth model with random se- 
quential dynamics that exhibits a transition from a 
smooth to a rough phase. In studying the model we have 



15 



shown that some properties may be understood directly 
from the scaling behavior of DP. However, other proper- 
ties and critical exponents appear non-trivial in the sense 
that they seem not to be related to DP quantities in a 
simple manner. Furthermore, we have introduced novel 
exponents such as those relating to the magnetization 
like order parameters characterizing the SSB and the re- 
sponse of these order parameters to an ordering field. 

One can also think of the model in terms of a system 
of unidirectionally coupled DP processes. We have pro- 
posed a field theory which should describe the properties 
of this general class of systems. 

We are left with several open questions. Firstly can 
the values of the novel exponents be predicted? Secondly 
can the critical behavior be described by a conventional 
scaling picture? In our study we have indicated that the 
critical behavior may be quite complicated. This could 
be consistent with multicritical behavior foun d in a study 
|p2fl of the field theory proposed in Sec. VII. 



We have also made a subtle connection between the 
present random sequential models and models similar in 
character to the parallel update PNG models such as 
that studied in H], It would be instructive to explore 
this point further. 

It would also be of interest to study in more detail the 
model (HH) that exhibits the breaking of a continuous 
symmetry. In particular the order parameter (|3^) has 
not been fully investigated. 

The growth models and their scaling behavior were 
investigated in one dimension. However, it is straight- 
forward to define the models in higher dimensions where 
similar scaling behavior is expected to hold. 

Let us finally remark that it would be very interest- 
ing to search for experimental realizations of the growth 
processes discussed in this paper, in particular because of 
their relation to DP. As pointed out by Grassberger Q , 
the large body of theoretical work on DP seems to be un- 
balanced by the fact that there are no experiments where 
DP exponents have actually been measured. It is not yet 
clear whether this is due to a lack of such experiments 
or to an oversimplification of nature in DP models. The 
growth models, however, suggest another category of ex- 
periments where DP exponents may be identified, namely 
absorption-desorption processes where the evaporation of 
atoms from the middle of completed layers is highly sup- 
pressed. 
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APPENDIX A: AN EXACTLY SOLVABLE CASE 

It is instructive to consider a case where the steady 
state can be calculated exactly. This will allow us to ver- 
ify that the interface is indeed rough in the moving phase 
and that if we allow some rate for desorption from the 
middle of a plateau, the interface is rough independent 
of whether the velocity is positive, negative or zero. 

The transition rates we consider are that of the RSOS 
model presented in the introduction with the addition of 
a process with rate p. 
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The process with rate p translates to desorption from the 
middle of a plateau when the model is translated back 
into a growth model via (0). We shall using a technique 
similar to that employed recently in jli] (see also p3[) to 
show that if 



P = 1 - .9/2 



(A2) 



the steady state probabilities P of a configuration may 
be written in a factorized form 



P(M) = Z N L 2 



1 o-M 



(A3) 



where M is the number of positive charges in the configu- 
ration and Zn is a normalization. One also has the con- 
straint that the only allowed configurations have equal 
numbers of positive and negative charges due to the fact 
that dynamics conserves the global charge. Taking this 
into account the normalization Zn is given by 



N/2 

£ 

M=Q 



-M 



M\ M\ (N — 2M)\ 



(A4) 



In order to 



prove Eq. (|A3j) let us define variables 

6-1-0, &++, b-i . . . where, for example, 6+o is the number of 

bonds where there is a positive particle at site i and 

a hole at site i + 1. Due to the fact that the global charge 

is zero we have b + o + b ++ + = feo+ + b ++ + b |_ = 

6_o + 6-- + 6-+ = &o- + b-- + 6_| = M which leads to 



b+o + b - + 26+_ 



"0+ 



6_o +26__ 



(A5) 



If Eq. ( |A3|) is to hold in the steady state the following 
equation for the balance of probability must hold for any 
configuration: 
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(6+o + gb 0+ + b - + gb- 
+26+_ + gb_+ + {g + p)b 00 )P(M) 
= (b 0+ + gb +Q + fe_ + gbo-)P(M) 

+(2 + g)b Q0 P(M + 1) + (, 9 6+_ + p6_+)P(M - 1) 



To understand Eq. (A6) note that the lhs gives the rate 
of loss of probability due to the transitions out of the 
configuration and the rhs gives the rate of gain of prob- 
ability due to the trans itio ns into th e co nfig urat ion. We 
now divide through by flAg) and use ( A5) in (]A6[) to yield 



In the case p = we see from (A2) that g — 2 and 
(A6) we are clearly in the moving phase of our original growth 
model. When g = 1 it is easy to check that we have de- 
tailed balance so that the interface is in eq uilibr ium and 
again we expect it to be rough. Equatio n (|A13 | ) verifies 
that in the case of non-zero p (and when (|A2|) holds) the 
interface is rough whether it be moving upwards (g > 1), 
downwards (g < 1) or not at all (g — 1). 



(g + P )b 00 - gb-+ = (1 + g/2)b Q0 + (2p - 2)6_+ , (A7) 

which is satisfied for arbitrary &oo, b |_ when ( |A2[ ) holds. 

In order to calculate the velocity and roughness for N 
large one notes that the sum for the normalization (A4) 
is dominated by M = N (1 — l/v2), therefore p, the 
steady state density of positive charges (and also that of 
negative charges), is 



(A8) 



p = i- i/V2 + e>(i/A0 



Here we define the velocity v as the steady state growth 
rate at an arbitrary site i. Let (ciCi+i) be the steady 
state expectation of finding a charge Cj at site i and a 
charge Cj+i at site i + 1, then 



g(0t - (+< Qi+i) 
-(0, -i+i) + g(-i Qi+i) 
+(3 5 /2-l)(0i !+ i) 
-2(+i - i+1 )+g{-i 



(A9) 



The form (A3) implies that correlation functions factorize 
to leading order in 1/N (e.g. (+j +i+i) — p 2 ) and one 
finds 

v = 2(g - l)p(l - 2p) 

+ (3.g/2 - 1)(1 - 2pf + (g - 2)p 2 + 0(1/N) 
~(2-V2)(ff-l). (A10) 

One can also calculate the roughness exactly. First we 
define the height at site i relative to site N as 



hi = c j 



(All) 



Clearly (hi 



0, therefore the width w is defined through 

N 

X 



1 N 



(A12) 



A moderate calculation then yields 
w 2 /N = p/3 + 0(1/N) ~ (1 - 1 



2)/3 



(A13) 



The result ( |A13|) implies that when ( A2 ) holds the inter- 
face is always rough and indeed the prefactor does not 
depend on g. 
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